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The singular value decomposition technique is used to reconstruct the electronic spectral weight 
function for a half-filled Hubbard model with on-site repulsion U = At from Quantum Monte Carlo 
data. A two-band structure for the single-particle excitation spectrum is found to persist as the 
lattice size exceeds the spin-spin correlation length. The observed bands are flat in the vicinity of the 
(0, tt), (tt, 0) points in the Brillouin zone, in accordance with experimental data for high-temperature 
superconducting compounds. 
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The study of the normal properties of the high- 
temperature superconducting compounds and the at- 
tempts to construct a microscopic theory of superconduc- 
tivity in these materials has been in large part dedicated 
to the investigation of the Fermi surface structure and its 
dependence on the strong electronic correlations and an- 
tiferromagnetic order. Angle-resolved photoemission ex- 
periments (ARPES) have recently provided a comprehen- 
sive mapping of the Fermi surface in some of the cuprates 
(YBa2Cu402 and Bi2Sr2CaCu208+x)|],[|,i. These 
results are particularly relevant to the high-temperature 
superconductivity scenario in which the phenomenon is 
related to the presence of van Hove singularities in the 
density of states close to the Fermi surface Quan- 
tum Monte Carlo simulations have recently been used to 
compare these experimental results with numerical calcu- 
lations in models of strongly correlated electrons 
An essential part of these numerical results was the ex- 
traction of the spectral weight function (SWF) for single- 
particle excitations from the Quantum Monte Carlo data 
calculated for imaginary time. 

The SWF is obtained as a solution of the following 
inverse problem : 

G{k,T)^ dLo _^p A{k,u)- (0<r</3) (1) 

where G{k,T) is the Matsubara Green's function for lat- 
tice momentum k and imaginary-time separation r, ob- 
tained in the finite-temperature Monte Carlo simulations 
at temperature and A(k,uj) is the corresponding 

SWF. 

This inverse problem admits an infinite class of solu- 
tions for A{k,uj) which will satisfy Eq.(|^) within small 
perturbations of the l.h.s., originating from the statistical 
noise on the simulation data. In recent works the one- 
electron SWF for the single-band 2D Hubbard model was 
obtained from Monte Carlo data using a maximum en- 
tropy approach The finite-size effects evident 
in the data itself, as well as in the reconstruction of the 



SWF, led the authors of Ref.]^] to conclude that the pseu- 
dogap in the single-particle spectrum observed in small 
clusters (up to 8^) is a lattice-size artifact which disap- 
pears as the size of the simulated system increases. They 
suggested therefore that for on-site repulsion weaker than 
the Hubbard bandwidth {8t) the antiferromagnetic cor- 
relation length would be comparable to the size of a suf- 
ficiently small lattice and would effectively create long- 
range antiferromagnetic order. Based on this conclusion 
the U = At coupling (where U is the on-site Coulomb 
repulsion and t is the nearest-neighbour hopping param- 
eter) was regarded as belonging to the weak-coupling 
regime, where the band structure is essentially similar 
to that of noninteracting electrons for all temperatures 
above zero, in accordance with the Mermin- Wagner the- 
orem. Following this argument, a gap originating from 
spin-density waves (SDW) will be both temperature and 
finite-size sensitive. The low-temperature simulations, as 
presented in for /? = 12 show indeed a stable gap even 
for relatively large lattice size. Such temperature depen- 
dence of the gap would be expected to correspond to a 
strong temperature dependence of the spin-spin correla- 
tions. However, as the numerical results of Refs.||ll[ and 
|]l2| show, the spin-spin correlations at t/ = 4i and 4^ lat- 
tice fall rapidly within 2 lattice spacings with only minor 
dependence on the temperature in the range /3 = 4 — 12. 
The SWF presented in Ref. j|] shows a clear two-peak 
structure for this value of U and lattice linear dimension 
larger than the correlation length. Hence the possibility 
that the band structure is due to the short-range antifer- 
romagnetic order (as was suggested in Ref.[||) demands 
an accurate mapping of the on-site repulsion regimes for 
coupHng U < 8t. 

A detailed study of the SWF as a solution of an inverse 
problem requires a quantitive criterion for the resolution 
limits of the reconstruction technique. Without such 
a criterion an existing gap might be overlooked in the 
reconstruction procedure. Moreover, the single-particle 
spectrum derived from the distribution of the function 
A{k,oj) along the Brillouin zone is based on the identifi- 
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cation of the peaks of this function. Thus the resolution 
of the two-peak structure, even if there is no clear evi- 
dence of the vanishing of this function between the peaks, 
is important. 

In our work we use a method based on the singular 
value decomposition (SVD) approach widely used in the 
field of inverse problems . This technique allows us to 
define the resolution limits quantitatively, thus providing 
an upper bound on the size of the gap with respect to 
the lattice size. 

The work presented is an investigation performed on 
lattice sizes ranging between 4^ and 12^. The simulation 
temperature was chosen as /3 = 5 and on-site repulsion 
as U = 4t. To avoid the sign problem the model was sim- 
ulated at half-filling. However, the SVD treatment can 
be extended to finite doping, as well as to other Green's 
functions. 

Eq.(^ can be rewritten in operator form as follows 



The solution of Eq.Q expanded in the functions {ufe} 
is given by: 



Gin) = {ICA){n). 



(2) 



K. is the integral operator of the r.h.s. of (|l|), with the 
variable r in the kernel discretised, acting on the func- 
tion A. This operator defines a transformation from the 
space of the SWF to the data space of which G{Ti) is 
a vector. This transformation may be regarded as a rect- 
angular matrix, /C, one of whose dimensions is infinite 
and the other equal to the dimension of the data space, 
namely to the number of data points rir- A conjugate 
operator K.* acts in the data space and is determined by 
the following equality of inner products: 



l^r.r,(/CA)(T,)G(r,) = J du;A{. 



c.)(/C*G)H, (3) 



with a metric tensor W which can account for the corre- 
lations in the data points. A common choice is to define 
this metric as the inverse of the covariance matrix . 
The operator /C/C* is represented by a symmetric square 
matrix with a positive set of n,- eigenvalues {af} and cor- 
responding set of orthonormal eigenvectors {vi}. These 
vectors form a basis in the data space. The operator IC*IC 
acting in the space to which the solution A(w) belongs 
has a set of rir orthonormal eigenfunctions {ui} with the 
same eigenvalues. If the set {cti}^JL^ is arranged in de- 
scending order each vector Vk and function Uk will have 
k — 1 nodes. These two sets {ui} and {vi} and the set 
{ai}"li form a singular system satisfying the equations: 



(4) 



The matrix representation of the operator /C/C* is ob- 
tained||l4| from the kernel of Eq.ffl) as: 



A{uj) = — < G,Vk > Uk{u}) 

k—1 



(6) 



(5) 



where <, > denotes the inner product in the data space. 
The smallest singular values will greatly amplify the con- 
tribution of the statistical noise to the expansion coeffi- 
cients in Eq.(^. We therefore introduce a cut-off into 
the expansion to exclude terms with ^ smaller than the 
statistical error in G(r) []l5| . 

The "resolution ratio" depends on the number of nodes 
in the function with the highest index kmax included 
in the truncated expansion [[l6| . Here prior information 
about the support of the solution is essential. Since we 
assume that the A{uj) is localised within a certain range 
of Lu, we confine the singular functions to this range. This 
is done by multiplying the kernel in Eq.(Q) by a pro- 
file function which is approximately 1 within the support 
and smoothly vanishes outside it . The rate of the de- 
crease of the singular values for the limited support grows 
as the profile function gets tighter, therefore, the trunca- 
tion in (|^) will have a smaller k^ax- However, since the 
support of the singular functions is limited, their nodes 
will be more closely spaced, thus giving better resolu- 
tion. To investigate the resolution limits, artificial data 
was generated by substituting an A{u!) consisting of two 
Gaussian peaks into Eq.(Q). The inverse problem was 
then solved for A{lu) with the expansion truncated 
at ^ exceeding a hypothetical noise level. Similar to 
the results of Ref. we found that for a fixed noise the 
resolution improves as the range of support decreases. 

The reconstruction of the SWF for k at the vicinity of 
k = (§, f ) is given in Fig.l. It can be seen that there are 
negative side-lobes. Just as a positive (5— function, im- 
aged over a limited Fourier bandwidth, turns into Airy 
rings which have negative side- lobes, we must expect 
a localised SWF, reconstructed over a limited singular 
function bandwidth, to have similar features. A regu- 
larization procedure which would restrict the solution to 
positive values only may also reduce the resolution of its 
reconstruction . We present here the results of the 8^- 
12^ lattices, since the existence of the gap in the smaller 
lattices is beyond controversy. Limiting the support of 
A{k, (jj) to the range [—8, 8], allowed us to identify a gap 
unseen in an infinite-support reconstruction. According 
to our statistical noise we truncated the expansion of 
A{k,uj) with the ratio ^ not exceeding 100 (for a 1% 
statistical error as estimated for 4000 Monte Carlo mea- 
surements) which allowed us to use 9 singular functions 
in the unlimited-support case and 7 in the reconstruction 
with the finite profile. 

Since the statistical noise determines the level of trun- 
cation of the singular-function expansion it is clear that 
at a particular truncation level, not all the data included 
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in the Ur points are of identical relevance. In accordance 
with the application of sampling theory in other inverse 
problems we identify as essential the nodes of the Wfc„„^ . 
This is a generalisation of the Nyquist rate in Fourier the- 
ory. For example, in the inversion of the Laplace trans- 
form |T7[ the optimal sampling was found to be expo- 
nential which corresponds to the nodes of Vk^^^ . Appli- 
cation of this sampling technique to the inverse problem 
(|l|) is performed as follows: we estimate kmax from the 
knowledge of the noise and the spectrum of the singular 
values. Then the data is evaluated at points nearest to 
the nodes of and the SVD reconstruction is done 

for these points only. This procedure leads to an obvious 
advantage from the computational point of view, since 
the evaluation of the Matsubara Green's functions is re- 
duced to ^ of the total points. Moreover, by separating 
the data points by a large spacing we reduce the corre- 
lations within the data and the metric on the data space 
can be taken as Euclidean, which also substantially re- 
duces the calculational effort. For the 8^ lattice we evalu- 
ated the SWF using 7 data points only and compared the 
results with that obtained from using 40 equally spaced 
points. For the latter case we accounted for the corre- 
lations using the covariance matrix as the metric in the 
data space. The results show remarkable agreement, thus 
fully justifying the sampling approach in this problem. 

In Fig. 2 we present the density of states calculated for 
the three biggest lattice sizes as the solution N{uj) of: 

/CO „-WT 

with G(0, t) being the zero-space-separation Matsubara 
Green's function. In this reconstruction the range of the 
support for N{uj) is [-12,12]. 

The gap decreases as the lattice size increases, in cor- 
respondence with the results of Vekic and White |8j , but 
it is still clear even in the 12^ lattice. In fact the dif- 
ference in the size and depth of the gap is very small as 
the lattice size increases from 8^ to 12^, while the lin- 
ear dimension of the lattice is considerably larger than 
the SDW correlation length. The existence of two clear 
peaks in the density of states (as well as in A{kFermi, ^)) 
suggests a two-band structure of the spectrum. 

The reliability of our results has been tested against 
the following criteria: the successful reconstruction of the 
Monte-Carlo data when the solution oi A{k,uj) is substi- 
tuted in the r.h.s. of Eq.(|l|); the successful reconstruction 
of the same data by: 

G=Y, < G,vk>vk; (8) 

k 

and, finally, by checking the first three moments of the 
SWF: 

fin= uj"A{k,uj)dLu (n = 0,l,2) (9) 



and the same moments for the density of states. We 
present in Table 1 the values of these moments for the 
density of states, showing excellent agreement with the 
analytic values 1,0,8 for the zeroth, first and second mo- 
ments respectively 1^. 



lattice 


Mo 


Ml 


M2 


8^ 


1.00 


0.00 


7.99 


10^ 


1.00 


0.00 


7.79 


12^ 


1.00 


0.00 


7.77 



Table 1. The first three moments of the density of 
states versus lattice size. 

Fig. 3 shows the single-particle spectrum derived from 
identifying E{k) as the location of a peak of A{k,uj), 
where the values of the momentum k scan over the Bril- 
louin zone. In the calculation of the spectrum we observe 
a significant difference between the 8^ and the 12^ re- 
sults. Both lattices exhibit a clear two-band structure of 
the single-particle spectrum. However, while the former 
can be fitted well with the mean-field SDW formula: 

Eik)MF = ^J{eikr + A^) (10) 

where e{k) — —2t{cos k^ + cos ky) and A = 0.64i, the 12^ 
lattice shows a bigger gap at the corners of the nonin- 
teracting Fermi surface than at the point (^, ^). Indeed 
the 12^ lattice simulation gives a smaller gap at this point 
than the 8^ and the lO'^ (see Fig.l). However, at the the 
points (tt, 0) the gap remains persistently clear for all the 
lattices. This result can be interpreted as suggesting a 
gap due to an anisotropic pairing channel Q| (although 
not necessarily of d^^-yi symmetry). Higher resolution 
calculations, based on data with better statistics will be 
required to examine this possibility. We also observe the 
feature of a flat band near k — (tt, 0) similar to the re- 
sults of Dagotto et al. [|j and the experimental data from 
ARPES @. This flatness may lead to a van Hove-like 
singularity in the density of states. 

The SVD method brings to this inversion problem a 
more rigorous and controlled approach than the maxi- 
mum entropy technique. Contrary to previous findings 
the results of our calculations suggest that U = 4i is in a 
regime in which the spectrum of quasiparticle (damped) 
excitations has a two-band structure, in spite of the on- 
site repulsion being sufficiently weak to allow substantial 
double occupancy. There is some support for the rel- 
evance oi U = At regime to high Tc superconductivity 
from a recent work by Beenen and Edwards [l8| . 
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P. Kornilovitch for stimulating discussions. This work 
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Figure captions 

Fig.la: A(7r/2,7r/2,w) for a 8^ lattice. Solid line - SVD 
solution with infinite support, dotted line - SWF support 

is limited in the range [-8,8]. 

Fig.lb: A(27r/5,37r/5,w) for a 10^ lattice. Solid line - 
SVD solution with infinite support, dotted line - SWF 
support is limited in the range [-8,8]. 
Fig.lc: A(7r/2,7r/2,w) for a 12^ lattice. Solid line - SVD 
solution with infinite support, dotted line - SWF support 
is limited in the range [-7,7]. 
Fig. 2: Density of states as a function ofu). 
Fig. 3a: 8^ lattice. Spectrum of single-particle excita- 
tions. Lower band - stars; upper band - circles. Dot- 
ted line - e(fc) ~ ~2t{coskx + cos ky); solid line - 

E{k)MF = e{kY + A^. Points in the Brillouin zone: 
r - (0, 0); M - (f , f ); M - (tt, tt); X - (tt, 0). 
Fig. 3b: Same as (a) for a 12^ lattice. 
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